Periodic ripples in suspended graphene 
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We study the mechanism of wrinkling of suspended graphene, by means of atomistic simulations. 
We argue that the structural instability under edge compression is the essential physical reason for 
the formation of periodic ripples in graphene. The ripple wavelength and out-of-plane amplitude 
are found to obey 1/4-power scaling laws with respect to edge compression. Our results also show 
that parallel displacement of the clamped boundaries can induce periodic ripples, with oscillation 
amplitude roughly proportional to the 1/4 power of edge displacement. The results are fundamental 
to graphene's applications in electronics. 



Graphene's unique electronic properties makes it ideal 
candidate for integrated circuits component. Theoreti- 
cally, conduction in a perfectly fiat graphene sheet can be 
ballistic, despite what observed experimentally is quite 
different. This is because that ripples are manifested^^ 
by developing a band gap, introducing additional effec- 
tive magnetic fields.^ Hence, understanding the rippling 
mechanisms is crucial for applications of graphene in na- 
noelectronics. Recently, it is reported that the ripple 
structure can be controlled by thermal treatment.!^ This 
brings out a straightforward way to the band gap en- 
gineering of graphene.^ Some relevant theoretical works 
have attempted to study this rippling using molecular dy- 
namics (MD), simulating suspended graphene under ax- 
ial compression at different temperatures.^ Graphene rip- 
pling by thermal treatment is known to be related to its 
negative thermal expansion coefficient (TEC),EIH1 and its 
membrane-natured mechanical properties.^ However, 
the intrinsic mechanism responsible for the periodic rip- 
pling of suspended graphene is not fully understood. 

From a mechanical point of view, there are in general 
two ways to induce periodic ripples in a suspended thin 
film: 1 Stretching in the axial direction which is perpen- 
dicular to the fixed boundaries, or, 2 Compression of the 
fixed edges in the lateral direction. Here we focus on 
the later case (so-called edge contraction), in view of the 
experimental observation of Bao et al.^ which shows bi- 
axial compression of graphene after annealing.^ We show 
the origin of edge compression due to heat treatment in 
Figjl](a): During heating, a difference in thermal defor- 
mation is created between the suspended graphene and 
the substrates, due to their different TEC. The graphene 
is stretched by the friction force given by the expand- 
ing substrate during the heating process. When, on the 
contrary, the system is cooling down, compressive force 
will be applied to the fixed boundaries from the interface. 
The deformation of the graphene boundaries becomes ir- 
reversible, because graphene exhibits high structure in- 
stability under in-plane compression.!^ 

In our MD simulation, we start with a graphene sheet 
suspended between two parallel supports, with an in- 
plane compressive strain imposed on its two edges fixed 
on substrate. The equations of motion are integrated by 
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FIG. 1: (Color online) (a) Schematic of edge compression 
of suspended graphene after heating and cooling, (b) Topo- 
graphic diagram of edge-compression-induced periodic ripples 
in a graphene sheet {wxL — ISnm x 5nm). Color scale shows 
the height profile (position along z axis). 



the Verlet algorithm with a time step of 1 fs. The system 
consists of 1000 — 20000 atoms. A Nose-Hoover ther- 
mostat is used to help the system to reach equilibrium 
at 300 K in the first 500000 simulation steps. We let 
the system progressively reach further equilibrium in the 
next 500000 steps without manual thermal control. Fur- 
ther details about the simulation techniques can be found 
elsewhere.!^ We note that the compressive strain is ap- 
plied on the fixed edges along the direction parallel to 
the trench, in order to generate ripples in experimentally 
observed direction. This makes the rippling orientation 
studied in this work perpendicular to that simulated in 
the relevant work of Abedpour et al.J We also note that 
the here-simulated ripples are different from the intrinsic 
ones due to thermal fluctuation, which were found to dis- 
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tribute randomly over the surface with height variation 
down to the atomic level 

The total interatomic potential involves many-body 
terms as a collection of that of individual bonds, 
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where and ip^ are the interatomic repulsion and at- 
traction terms between valence electrons, respectively, for 
bound atoms. The long-range interactions are included 
by adding a parameterized Lennard- Jones 12-6 po- 



tential term. 



tor 



is a single-bond torsion term. The 



bond order function hij includes the many body effects. 
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where h^~^ depends on the atomic distance and bond 
angle, hff represents the influence of bond conjugation. 
h^^^ is a dihedral- angle term for double bonds. This po- 
tential is an extension of the second generation of reactive 
empirical bond-order (REBO) modeLf^ The derivatives 
of —U^ (force components) are analytically computed. 
This semi-empirical approach has recently been used in 
many simulation works on the structural properties of 
carbon nanotube (CNTs) and graphene (e.g. Refs il2|17] - 
[2Q|) . It has also been used in one of our previous studies 
for investigating the nonlinear elasticity of CNTs,!^ in 
which we obtained quantitative agreement between the 
AIREBO-calculated Young's modulus of CNTs and that 
from ab-initio calculations. Compared to first-principle 
methods, an important feature of this empirical potential 
is its ability to deal with large systems. This is particu- 
larly important in case that the number of atoms cannot 
be reduced by using periodic condition. 

Fig|2]a shows an example of how the system attains 
thermodynamic equilibrium during the simulation, de- 
picting temperature and energy variations. We observe 
significant thermal fluctuations at the beginning of sim- 
ulation. During the thermal equilibrium phase, the tem- 
perature fluctuation is reduced from 150K to about 20K, 
which corresponds to the inter-atomic potential energy 
variation from 20 to 2.5 meV/atom. The simulation re- 
sult is taken as the average atomic configuration over 
the production phase. The resulting ripples are in well- 
ordered periodic wave shape (e.g. Fig TJ)), however, ran- 
dom disorder due to intrinsic thermal perturbation can 
be observed at any given instant during the simulation 
(Figj2|)). This is in agreement with the Monte-Carlo sim- 
ulation results of ReffMl 

The maximum number of atoms which can be simu- 
lated by our MD code is around 20000, which is already 
a large number in the atomic simulation world. How- 
ever, the graphene size used in experiments can be up 



FIG. 2: (Color online) (a) Temperature and potential energy 
variations during a simulation for the graphene shown in Figjl] 
(b). (b) Snapshot of the system at an arbitrary instant during 
the production phase. 



to 10^ times larger. It is therefore important to under- 
stand the influence of graphene size on the ripple struc- 
ture. In Fig Sja) we show the atomic configuration of 
a graphene sheet submitted to edge compressive strain 
e = Aw/w = 0.05. An interesting phenomenon observed 
in this figure is that the waves tend to merge into each 
other when they propagate from the edge to the center, 
e.g. 6 wave undulations can be observed near the edge 
(solid line in Figjsja), while only 2 are left in the sheet 
middle (dashed line). Also, it is found that the wave am- 
plitude in the middle is larger than that near the edges. 

To show the size effect, we plot in FigJsJ^) the height 
profile of ripple shape at the middle of three graphene 
sheets with different length L. We can see that, for a 
given edge contraction, the wavelength A and amplitude 
A increase with L. As a consequence, the ripples wave 
become less dense in the graphene middle and the wave 
form is changed. We find that the oscillation amplitude 
becomes about 4 times larger for the sheet 5 times longer. 
According to an analysis using the Foppl von Karman 
equations, the wavelength and the out-of-plane dis- 
placement in the region far away from the clamped 
boundaries should be roughly proportional to a/Z- In 
our simulation we observed that A and A both increase 
with L. The exactly linear dependence on ^/L is however 
not clearly shown, due to the fact that our graphene is 
not large enough. To show the full length dependence 
of the ripples, large-scale modeling approaches such as 
finite-element simulations are required for graphene size 
consistent with experiments, as pointed out by Shenoy et 
al.}^ 

Figjlja) shows the height profile at the middle cross 
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L= 24 nm 
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FIG. 3: (Color online) (a) Graphene (wxL — 12nm x 24nm) 
submitted to compressive edge strain e = Aw/w = 0.05. The 
solid lines show the shape of height profile, (b) Height profile 
of the middle cross section of graphene with different L. 



section of a graphene under compressive edge strains 
ranging from 0.01 to 0.1. It can be observed that more 
waves appear with larger oscillation amplitudes A when 
edge strain increases. Ref. 10 suggests that and 
should roughly hold a linear relationship with the lon- 
gitudinal strain 7. Since in our simulations the ripples 
are induced by a lateral deformation e , for establishing 
the correlation between the lateral strain e and A, we 
have used e ~ 7^; from the original definition of the Pois- 
son's ratio v. Hence, according to Ref.TO", the dependence 
of the ripple structure on the compressive edge strain e 
should be governed by 
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where t is the thickness. The Eqjs] suggests A 
and A ex e^/^. 

These linear dependences are clearly shown when we 
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FIG. 4: (Color online) (a) Height profile of the middle cross 
section of a graphene {w x L = 19Anm x 14.6nm) submitted 
to different edge strain e. (b) Wavelength A and out-of-plane 
amplitude A versus e^^^. The symbols represent simulation 
results and the lines stand for the best-fitted curves. 



agreement about the value of t of a one-atom-thick car- 
bon sheet, in particular when the surface is curved and 
involves Van der Waals interactions by tt- stacking. I^^l^ 
The data fluctuation in FigQ^) is due to the small size of 
the simulated graphene, since the 1/4-power law is only 
valid for the region far away from the clamped bound- 
aries. In the case of thermal treatment such as that 
shown in Ref|5l neglecting graphene's resistance to in- 
plane compressive stress, the strain e could be approxi- 
mated by 
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where AT is the magnitude of heating or cooling tem- 
perature, as and ag are TEC of substrate and graphene, 
respectively. Neglecting the second order term, we obtain 



plot the values of A and A as functions of e^/^. Tak- e ~ AT {as — (^g)- The wrinkling mechanism explained 

above also suggests that ripples can appear by only cool- 
ing without heating. 

Besides the thermal effects, our simulation results 
demonstrate that periodic ripples can be induced by in- 
troducing parallel displacement of the fixed boundaries. 



ing t = 0.339nm (as graphite's Van der Waals interlayer 
spacing), we find that the slopes of the two best-fitted 
lines suggest that the value of v is approximately 0.145, 
which is close to the ab-initio-calcuidited value. How- 
ever, it is necessary to note that there is still no general 
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FIG. 5: (Color online) (a) Atomic configuration of a sus- 
pended graphene {wxL — 13nm x 5nm) submitted to 0.65nm 
of edge displacement D (indicated by the arrow), (b) Out-of- 
plane amplitude A versus D. The symbols represent simula- 
tion results and the curve stands for the best-fitted one. 



As shown in Figjs] (a), in such a case the ripples are 
not perpendicular to the trench, as those observed in ex- 
periments without thermal-treated samplesP^The ripple 
orientation seems to be not very sensitive to the magni- 
tude of displacement D. The fitting curve of out-of-plane 
amplitude in Figjs] (b) suggests that A roughly follows a 
1/4-power scaling law of D sls A ^ 0.34D^/^. This cor- 
relation can be explained by the fact that should be 
proportional to the longitudinal strain 7/^ while the ef- 
fective value of 7 should be proportional to the edge dis- 
placement D. This results implies that the wavelength A 
is controllable by adjusting the displacement, since Xx A 
should be roughly constant. 

In summary, we have investigated the effects of edge 
contraction and displacement on the structural insta- 
bility of suspended graphene. Our simulation results 
show that periodic ripples are manifested in graphene, 
when longitudinal compression is applied to the fixed 
boundaries by thermal effects. The 1/4-power law for- 
malisms are found to be valid for the wavelength and 
the out-of-plane amplitude as functions of edge strain. 
The wavenumber is found to be larger at the graphene 
edge than that in the center. Besides the thermal effects, 
we show that parallel displacement of the fixed edge can 
induce periodic ripples in suspended graphene. In such 
a case the oscillation amplitude roughly holds a linear 
relationship with 1/4 power of the displacement. These 
results are essential for understanding the experimentally 
observed ripples on suspended graphene. 



^ A.L. Vazquez de Parga, F. Calleja, B. Borca, 
M.C.G. Passeggi, J.J. Hinarejos, F. Guinea, and R. Mi- 
randa, Phys. Rev. Lett. 100, 056807 (2008). 

^ A.T. NDiaye, S. Bleikamp, P.J. Feibelman, and T. Michely, 
Phys. Rev. Lett. 97, 215501 (2006). 

^ C. Lui, L. Liu, K. Mak, G. Flynn, and T. Heinz, Nature 
462, 339 (2009). 

^ A. Castro Neto, F. Guinea, N. Peres, K. Novoselov, and 
A. Geim, Rev. Mod. Phys. 81, 109 (2009). 

^ W. Bao, F. Miao, Z. Chen, H. Zhang, W. Jang, C. Dames, 
and C. Lau, Nature Nano. 4, 562 (2009). 

^ R. Miranda and A. 1. Vazquez de Parga, Nature Nano. 4, 
549 (2009). 

^ N. Abedpour, R. Asgari, and M.R. Tabar, Phys. Rev. Lett. 

104, 196804 (2010). 
^ N. Mounet and N. Marzari, Phys. Rev. B 71, 205214 

(2005). 

^ C. Lee, X. Wei, J. Kysar, and J. Hone, Science 321, 385 
(2008). 

^° E. Cerda and L. Mahadevan, Phys. Rev. Lett. 90, 074302 
(2003). 

S. Timoshenko and J. Gere, Theory of elastic stability 
(McGraw-Hill, 1963). 

V.B. Shenoy, CD. Reddy, A. Ramasubramaniam, and 
Y.W. Zhang, Phys. Rev. Lett. 101, 245501 (2008). 
^2 Z. Wang, M. Devel, and B. Dulmet, Surf. Sci. 604, 495 



(2010). 

'^^ A. Fasolino, J. Los, and M. Katsnelson, Nature Mater. 6, 
858 (2007). 

S. Stuart, A. Tutein, and J. Harrison, J. Chem. Phys. 112, 
6472 (2000). 

D. Brenner, O. Shenderova, J. Harrison, S. Stuart, B. Ni, 
and S. Sinnott, J. Phys.: Condens. Matter 14, 783 (2002). 
B. Ni, S.B. Sinnott, P.T. Mikulski, and J.A. Harrison, 
Phys. Rev. Lett. 88, 205505 (2002). 

H. Zhao, K. Min, and N.R. Aluru, Nano Lett. 9, 3012 
(2009). 

Z. Wang and L. Philippe, Phys. Rev. Lett. 102, 215501 
(2009). 

Z. Wang, M. Devel, R. Langlet, and B. Dulmet, Phys. 
Rev. B 75, 205414 (2007). 

Z. Wang, M. Devel, B. Dulmet, and S. Stuart, Fuller. Nan- 
otub. Garb. Nanostr. 17, 1 (2009). 

K.N. Kudin, G.E. Scuseria, and B.L Yakobson, Phys. 
Rev. B 64, 235406 (2001). 

Y. Huang, J. Wu, and K.G. Hwang, Phys. Rev. B 74, 
245413 (2006). 

T. Vodenitcharova and L.G. Zhang, Phys. Rev. B 68, 
165401 (2003). 

G.-G. Ghen, W. Bao, J. Theiss, G. Dames, G. Lau, and 
S. Gronin, Nano Lett. 9, 4172 (2009). 



